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Abstract 

Mixed-Integer  Programming  has  traditionally  been  restricted  to  deterministic  models. 
Recent  research  has  opened  the  door  to  stochastic  optimization  models,  which  are 
typically  dynamic  in  nature.  This  project  lays  the  foundation  for  stochastic  dynamic 
mixed-integer  and  linear  programming  (SD-MIP).  This  project  has  produced  several  new 
ideas  in  connection  with  a)  convexification  of  two-stage  mixed-integer  sets  and  b) 
multi-stage  (including  two-stage)  stochastic  linear  programming.  Together  a)  and  b) 
provide  the  foundations  for  SD-MIP  problems.  From  new  concepts  and  algorithms  to 
applications  and  software,  this  project  has  made  significant  breakthroughs  in  all  aspects. 
This  report  provides  a  synopsis  of  both  theoretical  and  computational  results.  As  a 
preview,  we  mention  that  currently  available  deterministic  MIP  solvers,  as  powerful  as 
they  are  known  to  be,  are  unable  to  solve  SD-MIP  models  of  modest  size  within  an  hour 
of  computing.  In  contrast,  our  decomposition  approach  provides  provably  optimal 
solutions  within  the  hour  time-limit. 

Introduction 

This  project  deals  with  the  solution  of  (SD-MIP)  problems  which  arise  in  several 
application  domains  of  interest  to  the  Air  Force.  For  instance,  multi- aircraft  trajectory 
planning  problems  with  aircraft  dynamics  lead  to  such  SD-MIP  problems.  Our  project 
focuses  on  the  mathematical  structure  of,  and  algorithms  for  SD-MIP  problems.  We 
begin  this  report  by  first  presenting  the  general  mathematical  structure  of  such  models. 

In  stating  the  SD-MIP  model,  we  borrow  the  notion  of  a  scenario  tree,  commonly 
used  in  stochastic  programming  (SP).  This  construct  (i.e.,  a  scenario  tree)  provides  the 
probabilistic  structure  of  information  that  may  be  expected  to  evolve  over  the  planning 
horizon.  A  scenario  tree  AT  consists  of  nodes  (n)  and  with  each  node  one  associates  a 
time  index  t(n).  The  links  of  a  scenario  tree  represent  the  flow  of  information  from  one 
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scenario  node  (n)  to  a  child  node  (n+),  and  one  has  t(n)  <  t(n+),  and  the  conditional 
probability  p(n  +  |n)  is  denoted  Pn+. 

An  important  concept  in  our  formulation  is  the  description  of  state  variables,  which 
consists  of  two  components,  that  is,  ^  :=  (xy^,  includes:  G  an  endogenous 

part  of  the  state  vector  (or  simply  “endogenous  state”)  \  and  ^  denoting  an  outcome  of 
the  exogenous  random  variable  (or  simply  “exogenous  state”).  Here  the  term 
“endogenous”  refers  to  the  fact  that  the  algorithm  may  be  used  to  exercise  some  control 
on  the  trajectory  of  (x^i),  while  the  term  “exogenous”  refers  to  the  fact  that  the  algorithm 
does  not  exercise  any  control  on  the  outcomes_^.  For  example,  in  a  stochastic  inventory 
(or  reservoir)  control  model  x^j  will  represent  inventory  (or  water-level)  while  ^  will 
represent  demand  (or  load/precipitation  resp.).  In  situations  where  we  need  to  refer  to  an 
entire  arbitrary  sample  path  (outcome),  we  simply  denote  it  as  m.  Finally,  we  adopt 
following  notational  convention:  The  root  node  is  indexed  by  node  0,  and  an  initial 
vector  Xq  (as  well  as  the  information  Sq)  is  assumed  to  be  given  (fixed). 

Given  a  scenario  tree,  the  SD-MIP  model  may  be  stated  as  follows. 

Min{d{[ito  +  £'[/Io+(5o+)]  '■  Uq  G  Uq,  Xo+  =  cIq+  +  Ag+Xg  +  5g+itg  }  (0) 


where  {/i„}  are  defined  recursively  for  n  >  1  as 

f  ^[^n+(£n+)]  1 

^71  Min  I  S.t.  tLfi  G  t/jj(Xjj),  i  (1) 

lUn(Xn)  =  [Un  G  M”""  X  j 

and. 


■^71+  ~  *^77+  +  ^71+^71-  (2) 

This  is  a  very  general  setting  for  stochastic  dynamic  mixed-integer  programming  because 
it  not  only  allows  a  multi-stage  decision  process,  but  it  also  allows  both  continuous  and 
discrete  variables  Un(.Xn)  —  ^  which  evolve  over  time. 


In  our  applications  (e.g.  mission  planning  for  UAVs  and  other  aircrafts),  discrete  choices 
(e.g.  choosing  which  “way-points”  to  use)  will  be  made  at  the  start,  and  revised  at  an 
intermediate  point  of  the  mission.  In  total,  this  constitutes  two  points  in  time  where 
discrete  choices  (way -points,  targets  etc.)  are  made,  and  several  time  points  at  which 
continuous  choices  are  enacted  (namely,  speed,  altitude,  etc.)  The  cost  associated  with 
carrying  out  such  missions  is  unknown  at  the  start,  and  the  total  cost  is  typically  revealed 
only  at  the  end  of  the  mission.  Thus  decisions  (both  discrete  and  continuous)  are  made 
under  uncertainty.  We  should  emphasize  that  our  models  are  best  suited  for  a  collection 


'  As  with  ^  we  use  a  uniform  size  of  the  vector  x^j  so  as  to  keep  the  notation 
manageable. 
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of  aircrafts,  so  that  the  synergy  can  be  utilized  to  get  the  most  out  of  the  available 
resources.  In  this  sense,  we  are  interested  in  vector-valued  decisions  and  constraint 
spaces. 


This  report  outlines  the  main  results  produced  to  solve  problems  whose  models  may  be 
looked  upon  as  Stochastic  -  Dynamic  Mixed- Integer  Programming  (MIP)  problems. 
While  the  presence  of  discrete  choices  leads  to  MIP,  it  is  important  to  recognize  that  the 
discrete  and  continuous  decisions  must  work  hand-in-hand.  As  a  result,  we  will  study 
both  stochastic  mixed-integer  programming  (SMIP)  as  well  as  stochastic  linear 
programming  (SLP)  problems.  By  using  a  combination  of  ideas  from  cutting  plane 
theory  of  deterministic  MIP  (especially  disjunctive  programming  of  Balas  1979)  and 
nested  decomposition  (e.g.,  Birge  1985,  Pereira  and  Pinto  1991),  the  algorithms 
developed  in  this  study  lead  to  several  important  breakthroughs. 

a)  In  the  SMIP  thrust,  we  have  developed  the  first  general  purpose  algorithm  for 
general  stochastic  mixed- integer  programming  problems.  Here  the  term  general 
is  intended  to  convey  the  notion  that  randomness  in  the  model  is  allowed  in  all 
data  elements,  and  the  discrete  decisions  are  allowed  to  be  general  integers,  rather 
than  simply  binary.  These  capabilities  are  particularly  important  because 
decisions  like  the  size  of  the  payload  to  deliver  matters,  and  they  are  typically 
discrete.  In  addition,  uncertainty  in  data  is  allowed  to  appears  in  every  type  of 
data  element,  be  they  objective  function  coefficients,  constraint  right-hand  sides, 
or  even  within  the  constraint  coefficients  themselves.  Prior  studies  of  this  class  of 
problems  have  usually  been  restricted  to  special  cases  (e.g.  binary).  While  others 
have  attempted  such  generality  in  the  past,  there  have  not  been  any  practical 
algorithms  with  which  to  conduct  computational  studies  in  the  case  of  general 
SMIP  models  considered  in  this  project.  The  algorithm  which  we  have 
developed,  referred  to  as  the  Ancestral  Benders’  Cutting  plane  (ABC)  algorithm  is 
by  far  the  most  general  approach  developed  to  date. 

b)  As  part  of  this  project,  we  have  also  developed  tools  for  very  large  scale 
Stochastic  Linear  Programming  (SLP).  There  are  several  reasons  for  this.  First, 
SLP  models  continue  to  challenge  many  of  the  fastest  computers  to  date,  and 
many  applications  within  the  DoD  (e.g.  planning  airlift  operations)  lead  to  some 
of  the  largest  optimization  models  to  date.  Towards  this  end,  we  have  developed  a 
new  notion  of  optimality  that  is  especially  geared  towards  very  large  scale  SP 
models.  This  new  concept,  which  we  refer  to  as  statistical  optimality,  allows  one 
to  predict  the  quality  of  a  proposed  decision  on  the  basis  of  statistical  tests  which 
measure  the  likelihood  that  a  proposed  decision  is  within  a  certain  tolerance  level 
from  an  optimal  solution.  With  this  new  concept,  we  have  been  able  to  show  the 
validity  of  our  algorithms  for  certain  instances  which  have  been  identified  by  the 
Defense  Science  Board  Report  (2011)  as  one  of  the  major  challenges  in  trade- 
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studies  conducted  for  DoD  acquisitions.  We  also  observe  that  just  as  LP- 
technology  is  critical  in  the  solution  of  large  deterministic  MIPs,  the  ability  to 
solve  large  scale  SLP  problems  is  just  as  important  in  the  solution  of  SD-MIP. 

We  begin  this  report  with  a  presentation  of  the  ABC  algorithm  for  the  two-stage  case,  and 
then  show  how  such  an  algorithm  may  be  applied  to  the  kind  of  mission  planning 
examples  that  were  mentioned  earlier.  Following  the  summary  of  the  discrete  case,  we 
will  present  the  ideas  underlying  our  SLP  advances.  In  all  cases,  we  provided  a 
illustrative  computational  results,  although  our  publications,  and  dissemination  material 
include  much  more  extensive  computational  results. 

Ancestral  Benders’  Cutting-plane  (ABC)  Algorithm  for  Discrete 
Choices 

The  discussion  presented  in  this  section  is  based  on  material  that  is  available  in  Qi  and 
Sen  (2015),  although  the  notation  used  in  that  paper  is  more  in  keeping  with  two-stage 
models.  Using  the  notation  introduced  in  (0-2)  a  two-stage  SMIP  model  may  be  stated  as 
follows: 

Min  [d^UQ  -I-  E[hn(Sn)]  ■  Uq  G  Uq,  -I-  B^Uq  ,Uq  G  x  }  (3) 

where  the  functions  {/i„}  are  defined  for  n  >  1  as 

f  dyi  Uyi  3 

KiSn)^clXn+  Min  I  S.t.  Un  e  Un(Xn),  i  (4) 

lUn(Xn)  =  [Un  G  E"""  X  DnUn  <  -  C„X„]  j 

where  the  probability  >  0  to  be  used  in  (3)  is  assumed  to  be  given,  and  the  number  of 
outcomes  indexed  by  n  is  finite.  As  in  (0-2),  the  initial  state  Xq  and  correspondingly  Sq 
are  given.  However,  in  contrast  to  the  statement  in  (0-2),  we  have  dropped  the  notation 
“0+”  in  the  objective  function  in  favor  of  using  scenario  nodes  with  index  “0”  in  stage  1, 
and  n  =  1, ... ,  iV  in  stage  2.  One  more  observation  about  the  second  stage  problem  in  (4) 
is  that  there  is  no  future  cost  because  the  second  stage  is  the  end  of  the  horizon.  This 
simpler  statement  in  the  objective  function  allows  us  to  focus  on  the  main  difficulty  at 
hand,  which  are  the  mixed-integer  requirements  in  the  constraints  for  Uq  in  the  first-stage, 
and  Uji  in  the  second-stage.  However,  we  will  return  to  the  case  of  terminal  (third-stage) 
costs  in  the  sub-section  which  extends  the  ABC  algorithm. 
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The  algorithm  consists  of  two  major  elements:  a  first  stage  branch-and-bound  (B&B) 
algorithm  in  which  the  objective  function  is  approximated  in  a  sequential  manner.  This 
B&B  scheme  is  an  extension  of  Benders’  decomposition  to  cases  where  the  objective 
function  is  difficult  to  approximate,  as  is  the  case  in  SMIP  models.  The  next  piece  of  the 
ABC  algorithm  is  to  provide  set-convexification  (of  the  second-stage),  and  value- 
function  approximations  which  are  easily  managed  within  the  B&B  process.  These  two 
elements  are  described  first,  and  we  follow  that  up  with  how  one  might  extend  this  idea 
to  the  class  of  DoD  applications  which  motivated  this  study.  We  also  provide  a  summary 
of  our  computational  experience  with  SMIP  instances,  many  of  which  are  unsolvable  by 
one  of  the  most  advanced  commercial  code. 

First  Stage  Branch-and-Bound 

The  ABC  algorithm  starts  out  the  branch-and-bound  process  (B&B)  process  at  node  0  of 
the  tree,  with  the  objective  function  cIqUq  -I-  £'[/iji(Sji)]  replaced  by  a  value  function 
approximation  as  is  common  in  dynamic  programming;  namely,  the  approximation  has 
the  following  form. 

Min  {djuo  +  ho(iio)  :  Uq  G  Uq.Uq  G  x  }  (5) 

As  is  common  in  the  MIP  literature,  we  assume  that  inputs  to  the  instance  include  upper 
and  lower  bounds  for  all  integer  variables,  and  these  box-constraints  will  be  denoted  Qq. 

While  value  functions  /ig  of  MIP  problems  are  notoriously  complex  to  discover,  and  their 
use  within  (5)  is  even  more  complicated  because  of  the  presence  of  mixed-integer 
variables  in  (5).  The  value  functions  /ig  have  a  close  kinship  with  MIP  duality,  and 
Caroe  and  Tind  (1998)  were  the  first  to  propose  using  MIP  duality  to  create  these  value 
function  for  use  in  (5).  However,  MIP  value  functions  are  themselves  non-convex,  and 
discontinuous  in  general,  and  solving  (5)  directly  is  impractical  at  best. 

Our  approach  is  to  replace  /lo(iig)  and  the  set  Ug  G  x  by  creating  functional 
approximations  /ig  y  that  are  valid  over  a  family  of  subsets  {Qqj},  which  cover  all  integers 
in  Qq.  Hence  the  lower  bounding  approximations  will  replace  (5)  by  the  following  family 
of  problems,  one  for  each  Q^  j. 

Min  [d^Uo  +  /i^j(Ug)  :  Ug  G  I/g,  Ug  G  Qlj  }.  (6) 

The  functions  [h^j]  will  be  required  to  satisfy  certain  properties  as  specified  below.  To 
do  so,  recall  that  —  an  +  A^Xq  -I-  B^Uq.  Hence,  for  a  given  Ug  one  can  associate  an 
MIP  value  denoted  /iji(sji(iig)).  Hence  we  will  require  the  approximations  ho  y(iig)  to 
satisfy 
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(7) 


/lo,;(Wo)  ^  ^  PnhniSn(Uo))  V  Uq  G  Qoj  ■ 

n 

Because  of  the  above  minorizing  property,  (6)  will  provide  lower  bounds  on  (5).  In  the 
B&B  phase,  each  subset  Q^j  will  be  associated  with  a  first  stage  decision  Uqj,  and  a 
lower  bound  Vq  j  defined  as 

v^j  :=  d^ulj  +  Wj)  :=  Min  {d^Uo  +  /ijy(uo)  :  Uq  G  Uo.Uq  G  Q^j  }.  (8) 

Then,  the  B&B  process  of  the  ABC  algorithm  will  estimate  the  most  optimistic  lower 
bound  as 


:=  Miny  [d^UQj  +  ho  j(uQ  j)],  and  G  argmin  {d^Uoj  +  /loj(ujy)],  (9) 

If  the  solution  Uq  satisfies  the  mixed-integer  requirement,  then  we  update  the  incumbent 
and  the  upper  bound,  as  in  any  B&B  algorithm.  On  the  other  hand,  if  Uq  does  not  satisfy 
the  mixed-integer  requirements,  we  identify  which  of  the  remaining  boxes  need  further 
exploration,  and  depending  on  the  exploration  rule,  one  would  refine  one  of  the  boxes  in 
the  collection  {Qoj}-  We  should  bear  in  mind  however,  that  further  branching  should 


also  be  accompanied  by  obtaining  improved  approximations  of  the  collection  of  functions 


Improving  Second-Stage  Approximations 


One  of  the  main  insights  underlying  the  ABC  algorithm  is  value  functions  to  be  used  in 
the  first  stage  can  be  improved  by  using  parametric  cutting  planes  similar  to  the  ones  that 
have  been  developed  in  the  context  of  stochastic  combinatorial  optimization  algorithms 
as  in  Sen  and  Higle  (2005)  and  Sen  and  Sherali  (2006).  One  of  the  major  advantages  of 
using  parametric  cutting  planes  in  the  second  stage  is  that  they  are  designed  to  adapt  to 
adjust  to  changes  in  the  first-stage  decision  which  is  a  parameter  for  the  second-stage  in 
the  context  of  decomposition  algorithms.  As  a  result,  these  cuts  automatically  provide  a 
memory  of  past  visits  to  fractional  solutions.  In  addition,  the  strengthened  relaxations 
provide  stronger  value  function  approximations  (i.e.  Benders’  cuts)  which  are  used  to 
update  {h^  j}  as  the  iteration  counter  k  increases.  These  second-stage  approximations  are 
described  next. 


In  the  notation  used  below,  each  node  nwill  have  variables  ii„,  and  these  will  be 
influenced  by  decisions  made  by  the  parent  node,  designated  as  node  0  in  the  two  stage 
case.  Because  we  are  interested  in  generalizing  the  two-stage  algorithm  to  the  multi¬ 
stage  case,  it  will  be  convenient  to  refer  to  the  parent  node  of  n,  as  the  node  n  — .  Thus, 
in  the  following  development,  the  variables  Uq  and  those  referred  to  as  will  be  the 
same.  This  shift  in  notation  will  facilitate  our  development  for  dynamics  after  stage  2  in 
our  multi-stage  application. 
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Because  multi-term  disjunctive  sets  have  been  shown  to  provide  the  convex  hull  of 
mixed-integer  sets,  we  adopt  the  use  of  multi-term  disjunctions  to  obtain  cutting  planes  to 
convexify  the  second-stage  problem.  These  cuts  can  be  generated  in  one  of  two  ways:  a) 
by  generating  disjunctions  based  on  the  cutting  plane  tree  (CPT)  algorithm  (Chen  et  al 
2011)  which  was  one  of  the  fundamental  results  obtained  from  the  2011-12  portion  of  our 
AFOSR-funded  research,  or  b)  using  a  truncated  B&B  tree  where  the  second-stage 
approximation  for  scenario  node  n  is  solved  until  either  a  maximum  time  or  B&B  node 
limit  is  reached.  At  that  point  in  the  algorithmic  process  for  the  second  stage  will 
typically  have  a  fractional  solution  which  can  be  cut  away  by  using  a  parametric  cutting 
plane.  The  resulting  cuts  will  then  be  added  to  improve  the  approximation  used  in  the 
second- stage  constraints.  The  details  of  this  strategy  are  summarized  next. 

Let  Un_  be  given,  and  suppose  that  we  proceed  to  a  scenario  indexed  by  n.  We  then  find 
an  approximate  solution  li^,  and  if  this  solution  satisfies  the  mixed- integer  requirements, 
then,  no  additional  cutting  planes  will  be  added  for  the  current  iteration.  On  the  other 
hand,  if  some  integer  restrictions  are  violated,  then  we  generate  a  lifted  cut  in  the  space 
liji).  The  cut  will  delete  the  point  u^),  but  will  not  delete  any  mixed- integer 
feasible  pair  In  order  to  accomplish  this  task,  let  us  first  examine  the 

following  set. 

1/nj  —  {(.^n>^n)\^n^n  ~  ^n^n-i  H  01 

Un-  e  Qn-,j-^^n  ^  {Qnj)]  }, 


where  ^  [u^  (Qn.x)]-  Note  that  the  matrices  and  and  the  vector  from  (1)  are 
subsumed  within  matrices  and  respectively,  and  the  vector  .  Also,  Xq  is  fixed. 
Thus,  substituting  the  state  variable  using  the  linear  dynamics,  we  obtain  a 
representation  completely  in  terms  of  the  variables  (Un-,UyJ.  This  set  has  the  following 
representation  (with  an  appropriate  matrix  and  the  right-hand  side  h^). 

^n,j  —  {(Mn-’^n)\^n^n-  "f  —  ^n>  ^n-  ^  Qn-,j’^n  ^  (Cn.r)]  }■  (H) 


The  set  in  (11)  is  a  disjunctive  relaxation  of  the  mixed-integer  requirement,  and  can 
therefore  be  used  to  generate  a  valid  inequality  using  disjunctive  programming.  Suppose 
that  the  resulting  cut  is  denoted  -I-  Then,  assuming  is  invertible, 

and  recalling  that  x„_  =  Xq  (fixed),  we  obtain  the  cut  coefficients  in  the  space  of  state 
and  decision  variables  and  the  set  (10)  can  now  be  updated  as  with 

appropriately  updated  matrices  and  the  vector  which  include  the  new 

parametric  cut. 
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Using  disjunctive  programming  results  (Balas  1979),  and  outcomes  of  our  AFOSR- 
funded  research  (Chen  et  al  2011),  we  have  shown  that  for  the  two  stage  case,  this 
approach  provides  a  finitely  convergent  algorithm  as  stated  below. 


Proposition  1.  (Qi  and  Sen  2015)  Suppose  that  the  set  of  first  and  second  stage  decisions 
form  a  compact  set  of  integers.  Let  the  cuts  be  derived  using  multi-term  disjunctive  cuts 
as  in  Chen  et  al  (201 1),  and  the  updates  of  the  approximations  in  (10)  and  (1 1)  are  carried 
out  accordingly.  Finally,  suppose  that  the  first  stage  B&B  method  uses  (6)  and  (7)  as 
discussed  above.  Then  the  ABC  algorithm  terminates  with  an  optimal  solution  in  finitely 
many  iterations.  ■ 

Application  of  SMIP  Models  for  Mission  Planning 

Consider  a  mission  dealing  with  carrying  out  reconnaissance  and/or  addressing  threats 
over  possibly  hostile  territory.  There  are  three  stages  of  decision-making,  with  the  first 
two  involving  the  discrete  choices:  a)  which  target/threat/way -point  to  assign  to  a 
particular  vehicle  b)  should  we  reassign  this  vehicle  based  on  observed  conditions,  and  c) 
identify  the  safest  route  back  to  base.  The  first  decision  will  be  denoted  Uq,  which  will 
specify  both  a  way-point,  as  well  as  the  pay-load.  While  the  way-point  and  pay-load 
combinations  are  important,  they  are  not  unique.  In  other  words,  if  at  some  point  in  the 
future,  the  way-point  changes,  there  may  be  a  subset  of  locations  which  can  be  served 
with  the  pay-load  decision  chosen  at  the  beginning.  As  in  the  presentation  of  the  ABC 
algorithm,  it  will  be  more  convenient  to  consider  this  variable  as  being  denoted  u^-,  with 
n  —  representing  node  “0”  of  a  scenario  tree.  Once  this  decision  is  made,  and  the  mission 
begins,  the  aircrafts  receive  data,  which  allows  them  to  make  some  adjustments  (recourse 
decisions),  which  will  be  denoted  ii„  where,  n  is  a  label  associated  with  the  specific  kind 
of  data  that  is  observed.  After  the  decision  is  implemented,  the  aircraft  will  carry  out 
the  mission,  and  will  use  a  sequence  of  maneuvers  which  will  deliver  the  pay-load,  and 
bring  the  aircraft  back  to  home  base.  These  decisions  will  be  represented  as  Un+  The 
decisions  associated  with  returning  to  the  base  may  be  looked  up  as  linear/quadratic 
programming  problem,  whereas  the  first  two  choices  are  based  on  discrete  selections. 
This  last  aspect  (delivering  the  pay-load  and  returning  to  home -base)  is  difficult  to 
predict  at  the  start  of  the  mission,  but  will  depend  on  which  specific  node  of  the  scenario 
tree  actually  transpires.  Hence  decisions  in  this  phase  will  be  denoted  u^+.  In  our  setup, 
the  decisions  ii„_  and  are  mixed-integer,  whereas,  ii„+  is  continuous.  Thus,  while 
the  model  is  in  fact  a  three  stage  model,  only  the  first  two  stages  have  mixed-integer 
variables,  and  the  setup  of  the  ABC  algorithm  work  is  applicable  to  this  situation.  A 
summary  of  these  choices  are  provided  in  Table  1. 
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Because  such  stochastic  dynamic  mixed-integer  problems  (SD-MIP)  are  computationally 
challenging,  it  is  important  to  trade-off  accuracy  with  tractability.  We  do  so  by  using 
linear  dynamics,  although  the  parameters  are  all  allowed  to  be  governed  by  observations 
which  could  be  stochastic.  Solving  such  a  model  helps  both  planning  and  execution 
because  optimal  (risk  minimization)  strategies  can  be  studied  prior  to  the  mission. 

Table  1.  Layout  of  Discrete  and  Continuous  Decisions  in  Mission  Planning 


Stage 

States 

State  Transition 
Equations 

Decisions 

(Continuous/Integer) 

0 

Xq  :=  Locations  of 
aircrafts  in  the  fleet 

Xq  is  given 

Uq  ■■=  Target,  Pay-load 
combination. 

(Discrete  Choices) 

1 

Xji  :=  Locations  and  data 
(conditions)  prior  to 
adjustments 

TAjjXq  -f  BjiIIq 

itji  :=  Re-assessment/re- 
assignment  based  on 
conditions.  (Discrete 
Choices) 

2 

x„+  :=  Locations  and 
data  immediately  after 
the  attack 

^n+  ~  ^n+ 
-\-An+Xn 

-f  By^^Uyi 

iiji+  :=  Return  path  to  base 
(Continuous  Choices) 

It  is  also  important  to  note  that  the  stagewise  approach  of  ABC  is  also  consistent  with 
many  multi-stage  SLP  algorithms  such  as  Nested  Benders  Decomposition,  Stochastic 
Dual  Dynamic  Programming,  and  others.  Accordingly,  this  application  is  very  well 
suited  for  a  combination  of  SMIP  and  SLP  models. 

Computational  Results 

Computational  results  with  the  ABC  algorithm  has  demonstrated  that  deterministic  MIP 
solvers,  even  the  most  sophisticated  ones  (e.g.  CPLEX)  are  not  reliable  enough  to  solve 
SMIP  models.  The  paper  by  Qi  and  Sen  (2015)  presents  many  examples  where  CPLEX 
fails  after  an  hour  computing  on  a  desktop  machine.  However,  the  same  machine  is  able 
to  compute  provably  optimal  solutions  using  the  ABC  algorithm.  Some  of  the  larger 
examples  reported  by  Qi  and  Sen  (2015)  are  summarized  in  Table  2.  Here  an  instance  of 
the  type  (mx5),  (nx5)  pertains  to  a  problem  with  m  potential  server  locations,  each  of 
which  can  accommodate  5  different  types  of  servers,  which  will  serve  n  potential 
locations,  again  with  at  most  5  different  types  of  clients.  As  one  can  observe,  CPLEX 
12.6,  which  is  one  of  the  most  sophisticated  deterministic  MIP  solvers,  and  has  been 
under  development  for  the  past  25  years,  struggles  to  solve  instances  with  a  large  number 
of  scenarios.  Interestingly,  when  CPLEX  12.6  does  find  an  optimal  solution,  it  does  so 
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within  fewer  CPU  secs,  (s)  as  shown  for  the  examples  in  Table  2.  However,  as  the 
number  of  scenarios  increases,  the  deterministic  solver  is  unable  to  scale  up,  and  failed  on 
6  out  of  9  instances  shown  in  Table  2.  These  results  demonstrate  the  practical  value 
behind  the  approaches  produced  as  a  result  of  this  project. 

Table  2.  Computations  Comparing  ABC  with  CPLEX  12.6  with  SMIP  Instances 


Instance 

\j\r\ 

Size 

Optimal 

Value 

Variables 

Constraints 

ABC 

(s) 

CPLEX 
12.6  (s) 

(2x5),(15x5) 

50 

-457.24 

1602 

1601 

4.26 

>3600 

(2x5),(15x5) 

100 

-420.35 

3202 

3201 

1.64 

0.78 

(2x5),  (15x5) 

500 

-399.47 

16002 

16001 

51.28 

>3600 

(3x5),(15x5) 

50 

-436.42 

2403 

1651 

191.92 

9.73 

(3x5),(15x5) 

100 

-371.24 

4803 

3301 

19.56 

>3600 

(3x5),(15x5) 

500 

-426.54 

24003 

16501 

1069.92 

>3600 

(4x5),  (10x5) 

50 

-308.69 

2204 

1201 

3.6 

0.98 

(4x5),(10x5) 

100 

-350.08 

4404 

2401 

261.89 

>3600 

(4x5),(10x5) 

500 

-311.3 

22004 

12001 

745.62 

>3600 

Large-scale  Stochastic  Linear  Programming 


As  mentioned  in  the  previous  section,  both  2-stage  and  multi-stage  SMIP  models  can  be 
embedded  within  the  setting  of  algorithms  which  solve  2-stage  and  multi-stage  SLP 
problems.  As  a  result,  we  are  also  interested  in  solving  large-scale  SLP  models  which 
can  end  up  as  a  platform  for  Stochastic-Dynamic  MIPs.  This  section  is  devoted  to  two- 
stage  and  multi-stage  SLPs. 


Statistical  optimality  for  two  stage  algorithms 

The  two  stage  SP  model  with  recourse  may  be  stated  as  the  following  optimization 
problem. 

f*  =  Min  /(uo)  :=  d(uo)  +  E[h(uo,oj)].  (12a) 

The  function  h  is  referred  to  as  the  recourse  function  and  is  defined  as 

h(uQ,(jL))  ■■=  Min{d7ui|  Du^  =  h(a))  —  C(a))iio,iii  >  0],  (12b) 

where  o)  denotes  an  outcome  of  the  random  variable  m  defined  on  an  appropriate 
probability  space.  When  the  space  of  outcomes  is  very  large  (or  we  have  continuous 
random  variables),  sampling-based  algorithms  provide  a  tractable  approach  to  two  stage 
SP.  For  such  algorithms,  optimality  can  only  be  verified  in  a  statistical  sense.  There  are 
two  classes  of  statistical  optimality  tests:  a)  In-Sample  stopping  rules  (e.g.,  Higle  and  Sen 
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1999)  and  b)  Out-of-Sample  stopping  rules  (e.g.,  Higle  and  Sen  1996,  Mak  et  1999  and 
Bayraksan  and  Morton  2011).  The  former  apply  to  algorithms  which  are  classified  as 
“internal-sampling”  algorithms,  while  the  latter  apply  to  both  “internal”  as  well  as 
“external”  sampling  algorithms.  The  current  project  has  designed  a  new  concept  called 
the  Compromise  Solution  for  out-of-sample  stopping  (Sen  and  Liu  2015). 

In  order  to  give  the  reader  a  concrete  sense  of  progress,  we  present  recent  computational 
results  with  the  Sonet  Switched  Network  (SSN)  cited  in  the  Defense  Science  Board 
Report  (201 1)  for  the  potential  power  of  SP.  This  is  a  two-stage  SP,  and  in  1994,  we  had 
reported  major  improvements  in  network  robustness  by  using  solutions  which  were 
obtained  on  a  Sun  Sparc  II  machine  after  about  an  hour  of  computing  using  SD  (Sen  et  al 
1994).  Because  the  SSN  model  involves  a  sample  space  of  0(10^*^)  possible  scenarios,  we 
had  no  hope  of  ascertaining  whether  the  design  we  obtained  was  close  to  optimum.  From 
a  deterministic  point  of  view,  this  is  true  even  today.  Even  if  one  had  access  to  exascale 
computing  (0(10  )  flops),  there  is  no  hope  of  predicting  the  expected  number  of  lost 
packets  deterministically.  If  one  accepts  statistical  assessments  of  optimality,  then  there 
is  hope. 

SAA  Computations  for  SSN 

The  methods  discussed  below  (i.e.  SAA  and  SD)  use  the  above  rules  in  a  manner  suited 
within  the  context  of  their  specific  algorithmic  procedures.  SAA  is  an  acronym  for 
Sample  Average  Approximation  (Kleywegt  et  al  2002),  and  the  strategy  here  is  to 
replace  the  expectation  E[h(u,oj)]  in  the  objective  function  of  (12a)  by  a  sample 
average  approximation  created  by  sampling  a  batch  of  N  outcomes  Thus,  the 

expectation  is  replaced  by 

N 

H(u)  =  — ^  h(u,  mO-  (13) 

t=i 

In  the  machine  learning  literature,  this  function  is  referred  to  as  “empirical  risk”.  Using 
concentration  inequalities  and  large  deviations  theory,  suggests  a  sample  size  which 
ensures  an  e-optimum  with  high  probability.  Nevertheless,  as  the  authors  of  SAA 
suggest  that  such  sample  sizes  are  much  too  large  to  be  computationally  practical 
(Kleywegt  et  al  2002).  Instead  a  sequential  application  of  SAA  is  recommended 
(Linderoth  et  al  2006).  For  each  sample  sizeiV,  the  sequential  approach  creates  m  = 
replications  denoted  of  the  form  (13).  After  optimizing  these 

M  replications,  one  calculates  a  lower  bound  confidence  interval,  and  an  upper  bound 
confidence  interval  as  suggested  in  Mak  et  al  (1999).  If  the  pessimistic  gap  is  large,  then 


^  The  pessimistic  gap  is  the  difference  between  the  upper  end-point  of  the  upper  bound  confidence  interval 
and  the  lower  end-point  of  the  lower  bound  confidence  interval. 
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the  sample  size  N  is  increased;  otherwise,  the  method  stops  and  declares  the  solution  to 
be  close  to  optimum. 

Table  3:  Results  for  SSN  using  SAA  and  LHS  Sampling  on  a  Computing  Grid 


Sample  Size 
(iV) 

Lower  Bound  (L) 
Confidence  Interval 

Upper  Bound  (I/) 
Confidence  Interval 

Pessimistic  Gap 
Max  (U  -  L) 

50 

10.10  (-h/-0.81) 

11.380  (-h/-0.023) 

2.113 

100 

8.90(-h/-0.36) 

10.542(-h/-0.021) 

2.023 

500 

9.87(-h/-0.22) 

10.069(-h/-0.026) 

0.445 

1000 

9.83(-h/-0.29) 

9.996(-h/-0.025) 

0.445 

5000' 

9.84(-h/-0.10) 

9.913(-h/-0.022) 

0.195 

The  results  for  SSN  in  Linderoth  et  al  (2006)  are  reproduced  in  Table  3.  For  N  = 
5000  the  pessimistic  gap  is  less  than  0.5%.  Thus  with  very  high  probability  the  optimal 
value  of  SSN  is  in  the  suggested  confidence  intervals.  The  results  in  Table  2  were 
obtained  using  a  computational  grid  with  several  hundred  desktop  PCs,  although  no  more 
than  one  hundred  machines  were  in  operation  at  any  given  time.  Even  so,  each  SAA 
instance  of  SSN  in  the  final  row  of  Table  3  (with  N  =  5000)  required  about  30-45  (wall 
clock)  minutes.  To  the  best  of  our  knowledge,  the  authors  use  M  —  6  for  SSN. 


Two-stage  SD  computations  for  SSN 

The  statistical  tests  used  in  SD  include  both  in-sample  as  well  as  out-of  samples  rules. 
The  in-sample  rule  verifies  whether  a  large  proportion  of  resampled  instances  provide 
duality  gap  estimates  below  a  given  tolerance  (loose,  nominal,  or  tight).  In  addition,  for 
out-of-sample  tests  SD  chooses  the  number  of  replications  m  =  Because  the 

specifics  of  the  algorithms  (SAA  and  SD)  are  different,  the  tables  present  slightly 
different  information,  although  the  bounds  and  pessimistic  gap  appear  in  both. 


^  Runs  with  N  =  5000  took  30-45  mins,  of  wall  clock  time  per  replication  with  100  PCs  on  a  computing 
grid 
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Table  4:  SD  with  Common  Random  Numbers  (on  a  MacBook  Air  with 

CPLEX12.4) 


Stopping 

Tolerance 

(Relative) 

Sample  Size 
(Standard 
Dev) 

Lower 
Bound  (L) 
Confidence 

Interval 

Upper 
Bound  (I/) 
Confidence 

Interval 

Pessimistic 

Gap 

Max 

(U-L) 

CPU  sec.  per 
replication 

Loose  (0.01) 

1023.33 

(167.62) 

9.366 

(-t/-0.244) 

9.953 

(-t/-0.050) 

0.881 

32.73(6.97) 

Nominal 

(0.001) 

2353.43 

(343.33) 

9.764 

(-t/-0.120) 

9.928(-t/- 

0.050) 

0.334 

109.96(26.31) 

Tight 

(0.0001) 

3137.50 

(605.17) 

9.876 

(-t/-0.107) 

9.925 

(-t/-0.050) 

0.206 

189.79(74.57) 

The  main  observations  from  the  above  comparison  (Tables  3  and  4)  are  provided  below: 

1.  (Choosing  decisions).  As  reported  by  Linderoth  et  al  (2006),  the  6  solutions  for 
SSN  found  using  SAA  were  quite  disparate  even  though  these  experiments  were 
done  using  Latin-Hypercube  Sampling,  a  variance  reduction  tool  commonly  used 
for  discrete-event  simulation.  This  was  attributed  to  the  inherent  ill-conditioning 
of  SSN.  Nevertheless,  solution  variability  (due  to  sampling)  does  appear  in  many 
practical  instances,  and  decision-makers  need  guidance  on  which  choice  to  adopt. 

2.  (Computing  Platforms  v  Algorithms).  The  MacBook  Air  has  a  processor  speed  of 
about  1.8  GHz,  and  the  average  Pentium  IV  processor  from  2004/2005  ran  at 
about  2  -  2.2  GHz.  Based  on  the  similarity  in  processing  speeds,  and  following 
the  arguments  in  a  PCAST  report  of  (Holdren  et  al  2011),  Sen  and  Liu  (2015) 
argue  that  the  speed-up  reported  in  Table  4  is  on  par  with  Moore’s  law,  as  long  as 
the  grid  computing  platform  spent  less  than  87.2%  of  its  time  on  either  idling  or 
communicating.  Similar  conclusions  have  been  drawn  while  comparing  Batch 
and  Online  Stochastic  Learning  (Bottou  2003),  where  “Batch”  refers  to  empirical 
risk  minimization  as  in  SAA,  and  “Online”  is  similar  to  SD.  Similarly 
computations  reported  by  Nemirovski  et  al  (2009)  also  support  this  conclusion. 
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Out-of-sample  Test:  Compromise  Decision. 

is  that  the  lines  on  the  left  represent  estimated  objective  function  values  at  the  solution  to 
a  Compromise  Problem  (14)  which  uses  the  replications  to  recommend  a  nearly  optimal 
decision. 


Mm  d(u)  +^\\u-u^  11^)  .  (14) 

V  m  / 

Consider  models  which  have  the  structure  shown  in  (12),  and  assume  that  i)  (i/,  fi)  are 
compact  sets,  U  is  a  convex  set,  and  d(u)  is  a  convex  function,  ii)  the  second  stage 
satisfies  the  relatively  complete  recourse  property,  and  iii)  the  random  process  affects 
either  activities  (columns)  of  the  first  or  resources  (right  hand  side)  of  the  second  stage. 
Let  an  interior  sampling  algorithm  (such  as  SD)  run  M  times  with  each  seed  creating  i.i.d. 
sequences  of  outcomes  {mjf},  form  =  1,  ...,M,  where  M  is  large  enough  and  N  denotes 
the  smallest  sample  size  discovered  by  SD  among  all  M  runs..  Suppose  each  run 
produces  a  solution  as  well  as  an  approximate  value  function  h^(u)  such  that 

G  argmm{/”^(ii)  d(u)  +  h^(u)  u  E  U] 

M  \  (15) 

and  \E[h(u^,a))]  —  =  Op  \Q^j- 

where  the  right-hand  side"^  in  (15)  is  due  to  the  Central  Limit  Theorem.  Let  Fj^(u)  denote 
the  grand  mean  value  function  approximation  —  Y,m  used  in  (1 1).  Suppose  that 

G  argmin{  Pm  (u) \u  E  U  ].  We  say  that  the  statistical  optimality  principle  is  satisfied 
if 


Fm(u^)  -  V  <f*  <  +  V>  (16) 

where  f*  denotes  the  optimal  value  in  (12),  and  p  —  Op  This  is  the  principle  that 

leads  us  to  the  kind  of  computational  performance  that  is  reported  in  Table  4.  Further 
details  for  SLP  problems  are  provided  in  the  Sen  and  Liu  (2015). 


Multi-stage  Stochastic  Linear  Programming 


This  part  of  the  project  deals  with  successive  approximation  schemes  within  the  context 
of  multi-stage  stochastic  programming,  and  should  be  interpreted  as  a  bridge  between 
stochastic  and  dynamic  programming.  The  method  presented  in  Sen  and  Zhou  (2014)  is 
the  only  MSLP  algorithm  which  discovers  the  stochastic  process  “on  the  fly”,  and  in  this 
sense  it  is  the  multis-stage  analog  of  two-stage  Stochastic  Decomposition  (SD).  Among 
sampling-based  approaches  for  MSLP,  the  most  popular  ones  trace  back  to  the  work  on 
stochastic  dual  dynamic  programming  (SDDP)  of  Pereira  and  Pinto  (1991).  Algorithmic 

Let  Op(l)  denote  stochastic  boundedness.  Any  sequence  of  random  variables  {R„}  is  stochastically 
bounded  if  for  every  e  >  0  there  exists  C  >  0  such  that  Pr[|R„|  >  C]  <  e. 
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enhancements  for  this  class  of  methods  have  appeared  in  Infanger  and  Morton  (1996)  and 
more  recently  as  abridged  Benders’  decomposition  in  Donahue  and  Birge  (2006). 
Asymptotic  results  by  Philpott  and  Guan  (2008)  highlight  the  need  for  certain  analytical 
safeguards  for  asymptotic  convergence. 

One  of  our  main  goals  was  to  present  an  approach  that  generalizes  the  two-stage  SD 
algorithm  to  the  multi-stage  setting.  More  specifically,  the  MSD  algorithm  provides  a 
dynamic  extension  of  regularized  multi-stage  SLP.  As  with  the  two  stage  case,  the 
multistage  extension  provides  an  asymptotically  optimal  solution  (wpl),  even  without 
requiring  the  scenario  probabilities  as  input  into  the  algorithm.  While  the  approach  for 
sampled  approximations  share  some  of  the  recursive  features  of  ADP,  some 
approximation  tools  (e.g.  duality)  and  asymptotic  analysis  draw  upon  concepts  that  are 
central  to  SP  (e.g.  epigraphical  nesting).  In  this  sense,  MSD  provides  a  bridge  between 
SP  and  ADP. 

In  order  to  appreciate  the  above  goals,  it  is  important  to  place  the  two-stage  SD  algorithm 
in  the  context  of  other  sample -based  methods  for  two-stage  SP.  Among  the  earliest 
approaches  for  sample -based  algorithms  we  have  the  Stochastic  Approximation  method 
(SA,  Robbins  and  Monro  1951).  The  work  of  Ermoliev  and  Shor  (1968)  appears  to  be  the 
first  application  of  SA  to  SP,  and  as  demonstrated  by  the  work  of  Nemirovski  et  al 
(2009),  there  is  continuing  interest  in  using  SA  methods  for  SP  problems.  This  genre  of 
methods  creates  a  sequence  of  sampled  subgradients  which  are  used  with  certain  step  size 
rules  to  ensure  convergence  (wpl)  of  the  random  sequence  of  solutions  generated  by  SA. 
Thus  one  might  consider  an  SA  algorithm  to  be  a  sampled  version  of  the  deterministic 
subgradient  method,  just  as  two-stage  SD  is  a  sampled  version  of  Benders’ 
decomposition.  The  advantage  of  these  sample-based  schemes  is  that  they  are  able  to 
work  with  statistical  estimates,  and  for  large  scale  models  they  have  the  advantage  of  not 
having  to  calculate  a  subgradient  using  every  potential  outcome  of  the  random  variable. 
The  main  point  here  is  that  SA  and  SD  are  both  randomized  versions  of  their 
deterministic  counterparts. 

In  order  to  extend  algorithms  that  use  Benders’  cuts  (e.g.  Nested  Benders,  SDDP,  etc.)  to 
cases  where  a  simulation  of  the  randomness  is  available,  but  not  exact  distributional 
information,  a  variety  of  issues  were  addressed. 

•  {Sampling  and  Convergence)  In  current  sampling-based  methods  motivated  by 
Benders’  decomposition  the  essential  role  of  sampling  is  to  reduce  the  number  of 
nodes  traversed  by  the  algorithm  in  any  iteration.  However,  the  approximations 
themselves  are  formed  by  using  deterministic  Benders’  cuts.  This  requires  a 
probabilistic  description  of  uncertainty  to  be  characterized  before  an  SP  algorithm 
can  be  used.  Note  that  there  are  many  emerging  applications  in  which  stagewise 
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independence  may  not  be  justified.  For  instance,  in  Philpott  and  Matos  (2012)  the 
authors  try  to  include  a  Markov  Chain  within  an  SDDP  framework  (see  also  Higle 
and  Kempf  2011).  However,  the  MSD  framework  provides  a  more  natural  setting 
for  such  applications  because  the  setup  is  based  on  a  dynamic  systems  framework 
which  admits  Markov  chains  seamlessly.  The  question  then  arises  whether 
asymptotically  convergent  algorithms  can  be  designed  for  such  applications.  For 
MSLP,  this  open  question  will  be  resolved,  in  a  mathematical  sense,  by  the  MSD 
algorithm.  This  paper  will  establish  asymptotic  results,  although  computational 
realizations  will  require  further  attention  to  data  structures,  stopping  rules,  and 
other  advances. 

•  (Uniqueness)  Under  what  circumstances  will  the  algorithm  produce  a  unique 
solution? 

•  (Specializations  for  stochastic  structure)  Are  there  specializations  (e.g.  stagewise 
independence  or  autoregressive  structures)  of  the  general  MSLP  structure  that  will 
allow  a  streamlined  algorithm  compared  with  the  case  of  correlated  random 
variables? 

•  (Specializations  for  constraint  structure)  Suppose  that  some  of  the  matrices  have 
specialized  structures  (e.g.  networks).  Are  there  ways  to  specialize  the  MSD 
algorithm  in  a  way  that  can  exploit  such  special  structures?  Finally,  does  the 
MSD  algorithm  reduce  to  two-stage  SD  when  applied  to  two-stage  SLP? 

By  developing  the  MSD  algorithm,  our  project  has  addressed  all  of  the  questions  posed 
above.  The  MSD  algorithm  is  a  dynamic  version  of  Higle  and  Sen  (1994)  which  is  a 
regularized  approach  to  the  original  version  of  SD.  In  the  SP  literature,  such 
regularization  was  already  an  effective  tool  within  the  deterministic  approach  proposed  in 
Ruszczynski  (1986),  as  well  as  non-smooth  optimization  methods,  as  in  Kiwiel  (1990). 
Finally,  these  ideas  were  also  studied  more  recently  in  Oliviera  et  al  (2012).  However, 
there  has  been  very  little  attention  paid  to  multistage  extensions  using  regularized 
algorithms  (e.g.  see  comments  by  Shapiro  2011).  The  MSD  approach  is  the  first 
regularized  algorithm  for  MSLP.  Details  of  this  approach  are  available  in  Sen  and  Zhou 
(2014). 
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1.  Dissemination 

Stochastic  Programming  Software  (Freely  Distributed) 


•  Available  at  www.neos-server.org/neos/solvers/slp:sd/SMPS.html 

•  Available  through  Github  at:  https://github.com/imliuvifan/sd 

•  Have  created  interfaces  for  PySP  (Open  Source  Software  for  Stochastic 

Programming) 

Presentations 

(Featured  Presentation)  S.  Sen,  “Stochastic  Mixed-Integer  Programming,”  INFORMS 
Optimization  Conference,  Miami,  FL,  February  2012. 

S.  Sen,  “Multi-stage  Stochastic  Decomposition,”  Symposium  of  the  Mathematical 
Optimization  Society,  Berlin,  August  2012. 

S.  Sen  and  Y.  Liu,  “Optimization  Statistics  for  Stochastic  Linear  Programming” 
INFORMS  Computing  Society  Meeting,  Santa  Fe,  NM,  January  2013. 

S.  Sen  and  Y.  Qi,  “The  ABC  algorithm  for  two-stage  Stochastic  Mixed-Integer 
Programming,”  International  Conference  on  Stochastic  Programming,  Bergamo,  Italy, 
2013  (also  presented  at  the  INFORMS  Annual  Conference  in  Minneapolis) 

S.  Sen  “Specializations  of  Multi-stage  Stochastic  Decomposition,”  Applied  Probability 
Society  Conference,  Costa  Rica,  July  2013 

S.  Sen  (Plenary/ Advance  Tutorial)  “Conceptual  v.  Computational  Stochastic 
Programming;  Perspectives  as  SP  turns  60”  INFORMS  Computing  Society  Conference, 
Richmond,  VA,  January  2015. 

Lectures  at  Universities 


90  min.  Lecture  at  Winter  School  at  Tignes,  France  (2013) 

Lecture  at  Ohio  State  University  (2013) 

Lecture  at  UC  -  Davis  (2014) 

Offered  ISE  638  to  Ph.D.  students  at  USC  (9  students  enrolled) 

Students  Supported 

Y.  Liu,  S.  Atakan  and  A.  Ozkan.  Their  dissertations  are  on-going  at  this  point. 
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Publications 


[1]  B.  Chen,  S.  Ku9ukyavuz,  and  S.  Sen,  “A  Computational  Study  of  the  Cutting  Plane 
Tree  Algorithm  for  General  Mixed-Integer  Linear  Programs,”  Operations  Research 
Letters,  Vol.  40,  pp.  15-19,  2012. 

[2]  D.  Gade,  S.  Ku9ukyavuz,  and  S.  Sen,  “Decomposition  algorithms  with  parametric 
Gomory  cuts  for  two-stage  stochastic  integer  programs,”  Mathematical  Programming, 
vol.  144,  pp.  39-64,  2012. 

[3]  K.  Huang,  S.  Sen,  and  F.  Szidarovszky,  “Connections  among  Decision  Field  Theory 
models  of  Cognition,”  Journal  of  Mathematical  Psychology,  Vol.  56,  pp.  287  -  296, 
2012. 

[4]  Y.  Qi  and  S.  Sen,  “Ancestral  Benders'  Cuts  and  Multi-term  Disjunctions  for  Mixed- 
Integer  Recourse  Decisions  in  Stochastic  Programming,”  revised  and  resubmitted  to 
Mathematical  Programming. 

[5]  S.  Sen,  “Stochastic  Programming,"  Encyclopedia  of  Operations  Research  and 
Management  Science  (S.  Gass  and  M.  Fu,  eds.).  Springer,  pp.  1486  -  1497,  2013. 

[6]  S.  Sen  and  Y.  Liu,  “Mitigating  Uncertainty  via  Compromise  Decisions  in  Stochastic 
Linear  Programming”  revised  and  resubmitted  to  Operations  Research 

[7]  S.  Sen  and  Z.  Zhou,  “Multi-stage  Stochastic  Decomposition”  SIAM  Journal  on 
Optimization,  Vol.  24,  pp.  127-153,  2014. 

[8]  H.  Gangammanavar,  S.  Sen  and  V.  Zavala,  “Simulation  and  Stochastic  Optimization 
for  Sub-hourly  Economic  Dispatch  of  Wind  Energy”  accepted  for  publication  in  IEEE 
Transactions  on  Power  Systems. 
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